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Abstract 

Coulomb blockade phenomena and quantum fluctuations are studied in meso- 
scopic metallic tunnel junctions with high charging energies. If the resistance 
of the barriers is large compared to the quantum resistance, transport can 
be described by sequential tunneling. Here we study the influence of quan- 
tum fluctuations. They are important when the resistance is small or the 
temperature very low. A real-time approach is developed which allows the 
diagrammatic classification of "inelastic resonant tunneling" processes where 
different electrons tunnel coherently back and forth between the leads and the 
metallic island. With the help of a nonperturbative resummation technique 
we evaluate the spectral density which describes the charge excitations of the 
system. From it physical quantities of interest like current and average charge 
can be deduced. Our main conclusions are: An energy renormalization leads 
to a logarithmic temperature dependence of the renormalized system param- 
eters. A finite lifetime broadening can change the classical picture drastically. 
It gives rise to a strong flattening of the Coulomb oscillations for low resis- 
tances, but in the Coulomb blockade regime inelastic electron cotunneling 
persists. The temperature where these effects are important are accessible in 
experiments. 

PACS numbers: 73.40.Gk., 73.40Rw., 71.10.+X 
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I. INTRODUCTION 



Quantum transport through mesoscopic metallic islands coupled to large reservoirs has 
been the subject of many theoretical and experimental investigations in recent years |1] [|. 
The small size of these systems implies a strong Coulomb interaction which gives rise to 
a variety of single-electron phenomena. When the temperature T and the voltage eV are 
low compared to the charging energy Ec, tunneling can be suppressed by the Coulomb 
blockade. However, when the energy difference Ao between two adjacent charge states is 
comparable to the temperature or the bias voltage, a current can flow through the system. 
As a consequence, the differential conductance shows a peak structure as function of an 
external gate voltage V g (linear response) or the bias voltage V (nonlinear response). In 
the absence of a dissipative environement and provided that the resistance Rt of a single 
barrier is much higher than the quantum resistance Rk = ^, i.e. ao = 4^2^ <S 1, the 
smearing of these oscillations is dominated by the temperature. In this classical regime, 
transport through the system is achieved by sequences of uncorrelated tunneling processes, 
which can be described by lowest order perturbation theory in the coupling between the leads 
and the metallic island. The corresponding classical rates can be used to set up a master 
equation H-fl, from which the charge average and the current can be calculated. For very 
low temperatures, however, or when the coupling «o becomes larger, quantum fluctuations 
set in fl9|-|i"5f and the classical picture breaks down for two reasons. 



First, the resummation of the leading logarithmic terms in aoln(J-^r) leads to a renor- 
malization of the difference in the charging energies, A , and the dimensionless conductance 
ao. Provided that A < T, we find 

A = A ° ~ = a ° m 

" ' '2a In (|§0 ' a ° 1 + 2ao In (!§0 " 

Throughout this work, we set h = k = 1 and consider only the case of a wide junction 
with many transverse channels. The renormalization is important if T < exp[— l/(2a )] ) 
which is an experimentally relevant temperature if ao is not too small. 

Secondly, coherent processes where the electrons can tunnel an arbitrary number of times 
between the leads and the island (resonant tunneling) give rise to an energy and temperature 
dependent broadening of the charge states. Thus, one can overcome the Coulomb-blockade 
not only by thermal activations but also by quantum fluctuations due to higher-order pro- 
cesses. For A ^> T this phenomenon is known as inelastic cotunneling [0,0] and can 



be described in second order perturbation theory in a . Important results of the present 
work are that higher order tunneling processes describe the crossover from cotunneling to 
sequential tunneling and can also give rise to important corrections to the classical result in 
the regime Ao < T where the conductance reaches its maximum value. We have called such 
processes "inelastic resonant tunneling" since, for a large number of transverse channels, all 
the electrons taking part in a coherent higher order process will be different. This is in con- 
trast to the usual mechanism of resonant tunneling through a single level in a quantum dot 
where only one electron is involved in any coherent process. We will show in this work that 
finite life-time effects due to resonant tunneling in metallic systems are especially important 
when 7r«o ~ 1- For finite temperatures this can be realized for a sufficiently large value for 
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In systems with low tunneling barriers we expect that effects due to energy renormal- 
ization and finite life-times are simultaneously observable in a real experiment at realistic 
temperatures. Both are important since wa <C 1 is equivalent to 2a o ln(^0 <C 1 except 
for very low temperatures which are experimentally not accessible. 

Charge fluctuations in the equilibrium case (e.g. in the single electron box) have been 
studied before by many authors. In ref. |13| a systematic perturbation expansion has been 
performed up to second order in ao including all possible charge states. This is a good 
approximation if the parameter 2aoln (^f) (and consequently also na^) is small compared 
to unity, i.e. the energy renormalization and the finite broadening are accounted for in 



a perturbative way. In ref. |ll],[T^,|rj[] the leading logarithmic terms together with certain 



improvements for larger values of «o |14[ have been considered in the two charge state 
approximation. However finite broadening effects are neglected. These approaches are valid 
in the low temperature regime where nao <C 1 and lead to the same renormalization effects as 
given by Eq. ([!]). In the nonequilibrium case (SET-transistor), the crossover from sequential 
tunneling to inelastic cotunneling has been studied in fl8|-f2T| by introducing a finite and 
constant lifetime into the expression of electron cotunneling. The results are valid in the 
parameter regime 2aoln(^^) <C 1 and tiolq <C 1 where renormalization effects and the 
energy dependence of the finite life-time can be neglected. In this regime several experiments 



have confirmed the theoretical predictions [ 22, 23 



In the present work, we develop a systematic diagrammatic technique in real time to 
identify the processes of sequential tunneling, inelastic cotunneling and resonant tunneling. 
If we concentrate on two charge states, we can resum the corresponding diagrams analytically 
and obtain closed expressions for the density matrix and all Green's functions. The spectral 
density describing the charge excitations of the system contains an energy renormalization as 
well as a finite broadening. Both are retained, which is crucial for a conserving theory which 
obeys sum rules and current conservation. The starting point of our technique, the real 
time representation of the density matrix, is closely related to path integral representations 
describing dissipation p4] , |2"5|j or tunneling in metallic junctions |5|j26[1. The method allows 
us to perform a nonperturbative analysis in the coupling to the reservoirs while taking into 
account exactly the strong correlations due to the Coulomb interaction. Thus, usual Green's 
function techniques, either for equilibrium |27| , [28[| or nonequilibrium |29|,|3(| systems, can not 
be used here. The same problem arises in the context of local, strongly correlated Fermi sys- 
tems like the Kondo or Anderson model [0l| -|33|| . For these systems diagrammatic techniques 
very similiar to our ones have been used by Barnes for the equilibrium case starting from 
a slave-boson description. Another example is a work of Rammer |35j who developed the 
same graphical language within a density-matrix description of the dynamics of a particle 
coupled to a heat bath. Finally, the technique presented in this work can also be formulated 
very elegantly in terms of Liouville operators using projection operator techniques developed 
in Ref. p6| . These relationships as well as the generalization to other systems like quantum 
dots with arbitrary many-particle correlations, time- dependent Hamiltonians and influences 
of the electrodynamic environement will be the subject of forthcoming works. 
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II. HAMILTON! AN AND PHYSICAL QUANTITIES 



We consider a small metallic island coupled via tunneling barriers to two leads. In 
addition it is coupled capacitatively to an external gate voltage (SET-transistor, see Fig. ([[])). 
This system is described by the following Hamiltonian 

H = H L + H R + Hj + V + H T = H + H T . (2) 

Here 

Hr — ^ e kn a krn a krn , Hj = £l n c\ n Cl n (3) 

kn In 

describe the noninteracting electrons in the two leads r = L,R and on the island where n is 
the transverse channel index which includes the spin. The wave vectors k and 1 numerate the 
states of the electrons for fixed r and n (a subindex k rn or l n has been omitted for simplicity). 
The Coulomb interaction V obtained by straightforward electrostatic considerations |I|,|4]Jl0| 
is given by 



V = E {^ - |) 2 (4) 



where E c = ^ is the charging energy, q x = C L V L + C R V R + C g V g , C = C L + C R + C g and 
V S ,C S (s=L,R,g) are the voltages and capacitances of the circuit according to Fig. ([[]). Qj 
denotes the charge operator of the island and is given by (e < 0) 

Q r = e(Nj - N + ) (5) 

where Ni = J2i n c \ n c in is the particle number operator and —eN + the positive background 
charge of the island. 

The charge transfer processes due to tunneling are described by 

H T = E J2(n>lnCln + h.C.) (6) 
r=L,R kin 

where T^ n are the tunneling matrix elements. 

For convenience, one can also introduce an auxiliary set of discrete charge states \N) 
with iV = 0, ±1, ±2, . . . and write the Coulomb interaction and the tunneling part in the 
form (see e.g. Q[E^£8|) 

V(N)=E C (N-^) 2 (7) 

e 



H t= E EW«Ue- !$ + h.c.) (8) 

r=L,R kin 

provided that we impose the constraint eiV = Qj in calculating the trace of any physical 
quantity. Thus eiV can be interpreted as the charge on the island which is now separated 
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from all the other degrees of freedom although they are of course not independent due to 
the constraint. $ is the phase operator canonical conjugate to N, i.e. [$, N] = i and e ±J * 
changes the charge by ±1 (note that N can take all integer numbers from — oo to oo here, 
so that the phase operator is well defined and the charge transfer operator e ±J * is unitary). 

Since both representations are equivalent, it is a matter of taste which one is used. 
However, for the metallic case the energy spectrum of the island is dense and the total 
particle number Ni on the island is very large. In this case the constraint eN = Qj can be 
disregarded since any canonical average over island operators with fixed Nj (or equivalently 
fixed N due to the constraint) can be replaced by the corresponding grandcanonical average. 
Such an approximation is of course not valid if we consider the case of a discrete level 
spectrum in a small island like a quantum dot. In this case, the particle number on the dot 
can be very low and one has to use the representation (f|) and @. 

Without the constraint, the separation of the charge degree of freedom is very convenient 



and effectively has also been used in other treatments like the path integral formalism p6 
The technical advantage is that the lead and island electron field operators occur now bilinear 
in the part Hq of the Hamiltonian (Q) since the Coulomb interaction (0) contains only the 
charge degree of freedom. In section 3 we will see that this allows for a straighforward 
application of Wick's theorem to integrate out all degrees of freedom of the leads and the 
island. 

The two leads as well as the island are treated as large equilibrium reservoirs which 
are not affected significantly by coherent tunneling processes where only a small number of 
electrons are involved. Thus we describe the electrons in these three "reservoirs" initially by 
Fermi distribution functions fi(E) = l/(exp(j3(E — fa) + 1) (i=L,R,I) with fixed chemical 
potentials fa or equivalently by the grandcanonical density matrix 



o - _L 

Pres r?Q eX P 
res 



/3 £ ( H i- PkNi 

i=L,R,I 



(9) 



where (3 = 1/T and Ni are the number of electrons in the leads and the island, respectively. 
In the following, we will always set the chemical potential \xj of the island to zero. 

Furthermore, we assume that at the initial time to there are no correlations between 
the reservoirs and the discrete charge states. This means that the initial density matrix of 
the total system reads po — Pres* where P° is a diagonal operator with matrix elements 
PjiJ- = (N\P°\N) and Pjy = P/v(to) is an arbitrary initial distribution for the discrete charge 
states at time t . 

After switching on the tunneling part Ht adiabatically, the probability P/v(£) f° r a certain 
charge N on the island at time t can be calculated from 

P N (t) = (N\e- iH(t - to) Po e lH{t - to) \N). (10) 

The stationary distribution follows from 

P$ = lim P N (t) = hm Pv(0) (11) 

t — >oo to — >— oo 

and will turn out to be independent of the initial choice for P^. Note that Pj$ is not an 
equilibrium distribution if the chemical potentials Hl/b, of the leads are different. 
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The average charge of the island, which is an experimentally measurable quantity via 
the electrostatic potential of the island, can be calculated from 



N = Y,NP^. (12) 

N 

The current flowing through the barriers r = L/R is defined by I r (t) = e^(N r (t)) . 
After a straightforward calculation one obtains 

l r {t) = 2elm lT,nr((a{ rn c ln e-^)(t))\ (13) 

I kin ) 

and for the stationary current 

If = lim iJt) = lim I r (0). (14) 

t^oo to— »— oo 

At the end of section 4 (see (|56|) ) we will show that the stationary current can also be 
related to the real-time correlation functions 

C>(t,t') =-i( e -*) e *)) /)o (15) 
C<(t,t') =i(e^ t ' ) e- i ^) po . (16) 

Both are independent quantities since we do not assume equilibrium here. In the stationary 
limit, i.e. to ~~ > — oo, the correlation functions depend only on the relative time difference 
t = t — t' and we can define the Fourier transform 

C>{uj) = J dre i£Jr C > (r) (17) 

and analog for C < {uj). 

From a theoretical point of view, an interesting quantity is also the spectral density 

A{u) = ±[C<{u)-C>{u)\ (18) 

which describes the charge excitations of the system. At the end of section 4 we will see 
that within our approximations all quantities of interest like the probability P$, the current 
I^i R and the correlation functions C > ,C < can be related to the spectral density. Thus, 
its specific form, which is related to energy renormalization and finite life-time broadening 
effects, will be the central point of our analysis. 



III. DIAGRAMMATIC TECHNIQUE 

To start with, we consider the probability distribution ( |10|) and write it in the form 

P N (t) =Y,P°N>Tr res p° res (N'\TMe l ^ dTHT{T)o \N) (N\ T^e^o ^ t(t)o \N') (19) 

N' 
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where Tr res is the trace over all reservoirs, T^ + > (T(~*) are the time-ordering (anti-time- 
ordering) operators and Ht(t)q denotes the tunneling part of the Hamiltonian (Q) in the 
interaction representation with respect to Hq. Note that the Coulomb interaction V is 
included here in Hq and will be treated exactly in the following. Eq.flTUD has also been 
used as a starting point by other authors, e.g. by Feynman and Vernon |24| or Caldeira 



and Leggett |!9| (see also for a system coupled to a heat bath or by Eckern et al |26 
for a tunnel junction. After integrating out the reservoirs within a real-time path integral 
representation these authors obtained an effective action where the two propagators occuring 
in Eq.(|Tj|) are coupled to each other. Expanding the exponential of this action in terms of 
the tunneling, one can arrive at a graphical language in real-time space || . 

In this work we will use a different approach which is easily generalized to other systems 
such as quantum dots, Anderson and Kondo models, etc. |40| . The procedure is first to 
expand the propagators in Ht 

T^e ±l ^ dTHTiT)o = £ (±*r fdn I* dr 2 ... r 1 dr n T^{H T (n) • ■ ■ H T (r n )} (20) 

and, in a second step, insert the form (§) for Ht and apply Wick's theorem with respect 
to the reservoir field operators a^} n and Cj„ . This is possible since H is bilinear in these 
operators. As a consequence, the vertices corresponding to the charge transfer operators e ±% ^ 
become coupled by reservoir lines, either from the two leads or the island. This is indicated 
in a graphical language by solid lines for the leads and wiggly lines for the island as shown 
in Fig.£| There the upper (lower) horizontal line corresponds to the forward (backward) 
propagator and the vertices change the discrete charge states as indicated by the particle 
numbers associated with the horizontal lines. The two charge states at the left end of 
the diagram at time to are the same and correspond to the initial probability distribution 
Pjyv (if nothing is written there in a certain graph we imply this automatically). The two 
charge states at the right end at time t are chosen equal to N if we want to calculate Pjv(t)- 
Graphically we indicate the relationship of the forward and backward propagator at the right 
end point symbolically by |iV)(jV|. Furthermore we assign to each line a certain energy, the 
Coulomb energy V(N) to the horizontal lines and the energies e (E) to the lead (island) 
lines. The rules to translate a certain diagram in time space are the following ones: 

1. Assign a factor e~ tx i T i to each vertex where Tj is the time variable and Xj is the difference 
of the energies entering the vertex minus all energies leaving the vertex. 

2. To each loop formed by 2s reservoir lines we assign a factor 

= e ^z im (^ wr 2??2 (^ e 2 )... r[~(6 S) e s ) (21) 

l\...ls 

where e±, . . . , e s and Ei, . . . , E s are the energies associated to the lead and island 
electrons, respectively, ordered in the direction of the loop, Tj = L/R specifies the 
left or the right lead and each sign o^r\j = ± indicates whether the reservoir line with 
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energy ej,Ej runs to a lower (+) or a larger (-) time with respect to the closed time 
path formed by the two propagators. The quantities V are defined as 

rSTfe £)=E TSTSTSie - el n )5(E ~ Zn)f?{e)f?W (22) 
k 

where = f\ (i=L,R,I) is the Fermi distribution function and ff — l — fi. 

3. The prefactor is given by (— i) M (— l) m (— l) c , where M is the total number of internal 
vertices, m the number of internal vertices on the backward propagator and c the 
number of crossings of reservoir lines. 

Finally one has to integrate over all time variables Tj from to to t and over all energy variables 
of the reservoir lines. 

For the current ( |TB| ) and the correlation functions (|TJ|[L1d|), the procedure is completely 



analog, the only difference is that additional external vertices have to be added. As shown 
in Fig.^| one has to introduce an external vertex with two reservoir lines at the right end of 
the diagram to calculate the current, and two external vertices to calculate the correlation 
functions. 

So far, the diagrammatic rules are formulated in time space which is not the most 
convenient one to calculate stationary transport properties. Therefore, using (|TT|) and (|i~4|) 
we set t = and to = — oo in every diagram for Pjv(t) and I r {t) and evaluate the time 
integrals analytically for a given fixed ordering of all time variables by using the identity 



f dr x f dr 2 ... [ dr M e~ ixiT1 e~ iX2T2 ■ ■ ■ e - iXM ™e VT1 

J — OO J T\ J TM — 1 



Xi+irj x 1 +x 2 + irj x x + x 2 + ■ ■ ■ + x M + iv 



(23) 



where e VTl (r/ — > + ) is a convergence factor which is related to the adiabatic turn on of the 
coupling part H?- The rules in energy space read 

1. For each auxiliary vertical line which does not cut any vertex and has always a vertex to 
its left, we assign a resolvent A ^ +ir) where AE is the difference of all energies crossing 
the vertical line from right to left minus all energies crossing it from left to right (see 
Fig.|). 

2. For each loop we assign the same factor as given by Eq.([2T|). 

3. The prefactor is given by (— l) m (— l) c , where m is the total number of internal vertices 
on the backward propagator and c the number of crossings of reservoir lines. Fur- 
thermore, if we are calculating the current, we have to multiply with an additional 
factor — ie and we have to assign a factor —1 if the external vertex to the right has an 
incoming lead line. 



For the Fourier transform ( |TTD of the correlation functions we write 



C>(cu) = / dr 



o 



e~ iuJT lim C > (0,r) + e^ r lim C>(r,0) 

£q — >— oo to — >—oc 



(24) 
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and again calculate all time integrals analytically by using ( p3[) and fixing the time ordering 
of all time variables including the integration variable r. Rule 1 has then to be supplemented 
by 

4. If an auxiliary vertical line lies between the two external vertices, we have to add ±uj 
to the energy difference AE when the vertex e ±l ^ lies to the left of the auxiliary line 
(see Fig.§). If we imagine a virtual line which connects the two external vertices from 
e~ % ^ to and assign an energy w to it, we can simulate the same by consequently 
applying rule 1 including this virtual line. 

The graphical representation of C > (u) and C < (u>) in energy space is given in Fig.|| 



IV. PHYSICAL PROCESSES 

Usually, tunneling processes within the SET-transistor are described by a classical master 
equation §§-§0 



P$(t)= E [P$(t)lN'N-P$\thNN'] (25) 

where ^n'n are classical transition rates calculated by the golden rule in second order per- 
turbation theory in Ht 



JN'N 



2tt^2 [ a r (Ajv)5jv,jv+i + at(^N-i)S N ',N-i • (26) 



Here Ajy = V(N + 1) — V(N) is the change of the Coulomb energy if we increase the excess 
particle number by one and 

a±(u) = J dE a±' T (u + E, E) (27) 

are the transition rates for tunneling in and out processes. Assuming constant tunneling 
matrix elements T rn = T£™ and neglecting the energy dependences of the density of states 
Pr ( e )> Pi(P) i n the leads and the island, we obtain the well-known expressions 

af(u) = {u- fi r ) o r n t (28) 

where 



Or. 



1 f? 

o = El^ \PrPi = ^~2^r (29) 



is the ratio of the quantum resistance Rk = \ to the tunneling resistance of barrier 
r, n+{oj) = l/(exp(P(u — fi r )) — 1) is the Bose function and n~ — 1 + nf. Furthermore, 
we introduce the notations ao = J2r a o > o ± {oj) = J2r 0^(00) , a r (uj) = a^{oj) + ol~(u) and 
a{u) = a + {oj) + a~{uj). 
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Physically, the master equation (|25|) describes sequences of uncorrelated lowest order 
processes (sequential tunneling) where each single process describes one electron entering or 
leaving the island. Retardation effects and higher-order correlated tunneling processes are 
neglected within this approach. We can write (p5|) in the form 

P { N\t)=Y. P Nht)lN>N (30) 
N' 

with 7tv'at = 7at'at — ^jv'JV J2n" Jnn"- We will now show that the exact kinetic equation has 
a similiar form 

*W) = E f dt' P N ,{t')t N , N {t' -t). (31) 

S denotes the sum of all possible correlated processes. Now non-Markovian effects are 
included. For the stationary solution Pjy this gives 

= y^P^/Sjy/jy , (32) 
N' 

where Eat/at = i f® 00 dr'EN'N(T r ). This equation can immediately be obtained by using our 
diagrammatic technique in energy space. The self-consistent equation for Pjft is shown 
in Fig.|7|. Here £ is defined graphically as the sum of all possible irreducible self-energy 
diagrams, which are defined such that any vertical line cutting them will cross some reservoir 
line (see Fig||). According to our diagrammatic rules 1,2 and 3 in energy space, the equation 
corresponding to Fig.^ reads 

Pn — P*n + Pn'^n'n — - (33) 

Multiplying with irj and performing the limit rj — > + , we arrive at ([32]) and can identify 
the sum of all possible correlated tunneling processes diagrammatically by the irreducible 
self-energy S. Furthermore, the solution of Eq.(p2[) is independent of the initial choice for 
Pjy which has dropped out of Eq.(^) in the limit rj — > + . In the same way, one can also set 
up a self-consistent equation in time space to get the time-dependent kinetic equation (ffTp. 

By moving the last vertex to the right of each diagram of £ up or down, one can easily 
show that 

E S ^v' = ° (34) 

N' 

which provides the possibility to write Eq.(^) also in the form 

= ^2 (Piv'^N'N — Pn^nn'J ■ (35) 

jvviv 

This has the form of a balance equation for each charge state \N). Furthermore, by changing 
the vertical position of all vertices and the direction of all reservoir lines one can easily see 
that Xat'at is purely imaginary, i.e. there exists a real solution of Q3"2]) and fl3"5|). 
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Calculating Eat'tv in first order in a (see Fig.|9|), we get 



s 7v+i,7v = 2lT ia (Ajv) , ^n-in = 27ria + (AAf_i) , (36) 



which, after insertion in Eq.(p2[), leads to the classical result ([30|) in the stationary limit. 

Before starting to include higher-order processes in £, we will introduce two important 
approximation which will simplify the following analysis considerably. 

First we assume that the number Z of transverse channels is very large. Since each loop 
contribution (|2lD is proportional to ZT S , this has the consequence that all loops with s = 1 
will dominate in each given perturbation order in T. Thus, we can restrict ourselves to 
the loops with two reservoir lines which give the contribution a^ ,T (e, E). Since the rest of 
the diagram depends only on uj = e — E, we can integrate af ,T (oj + E, E) over E and get 
the contribution af(uj) defined by fl2Tp8|) . Furthermore, we represent all loops with two 
reservoir lines from now on graphically by single solid lines with one energy variable uo (see 
FigjTJ. 

Secondly, if the charging energy Ec is large compared to the energy difference of the two 
charge states A = V(l) — V(0) and the bias voltage V — Vl — Vr, we can use the two 
charge state approximation. This means that we consider only the charge states N = 0,1. 

The classical approach breaks down for low temperatures or for large values for the 
dimensionless conductance ao- Specifically we will see in section 5 that the classical master 
equation is valid for a In (^^) <C 1. To go beyond this regime we consider now higher-order 
correlated tunneling processes in E. Fig .[11] shows an example in second order in «o which 



represents a part of the rate of inelastic electron cotunneling [16,17 when the leads r and r' 
are different. Here one electron enters the island from lead r, the system stays in a virtual 
intermediate charge state with N — 1 and finally another electron leaves the island via lead 
r' . However, the calculation of all second order diagrams for S i = — S o or S 10 = — £n 
is plagued by several irregular integrals which occur since the complete self-energy depends 
nonanalytically on a - Thus, we will directly perform a nonperturbative resummation of 
higher-order diagrams and discuss the second order result as a limiting case at the end. 
By this, we are also able to clarify certain singularities which are present in the usual 



expressions for inelastic electron cotunneling |16|,17j|. Furthermore, we are able to go beyond 



the cotunneling theory and can investigate the influence of resonant tunneling processes 
which will modify the classical result of sequential tunneling significantly if «o In (^0 ~ 1- 
A very illustrative example of a resonant tunneling process is shown in Fig.[l2|. Here the 
charge of the island is alternating between and 1 via an infinite number of intermediate 
virtual states. The electrons going back and forth between the leads and the island or 
tunneling from one lead to another via the island are all different since we have assumed a 
large number of transverse channels. This is in contrast to the usual mechanism of resonant 



tunneling where only one level of the island is involved. Of course, the diagram of Fig.|12 
is not the only one which is important to describe higher-order processes. In order to get 
a systematic criterion which diagrams might be relevant ones, we take a look at the states 
corresponding to the upper and lower horizontal lines before integrating out the reservoir 
degrees of freedom. Two such states which lie on the same vertical cut define a matrix 
element of the total density matrix (i.e. reservoirs plus charge states). Our criterion now is 
that we take into account only those nondiagonal matrix elements of the total density matrix, 
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which differ at most by one electron-hole pair excitation in the leads or (equivalently) in the 
island. Graphically this means that any vertical line will at most cut through two solid lines 
each representing a pair of one lead and one island line (see Fig.|13| for examples). The sum 
of all these diagrams can be represented as shown in Figs. |I~4] - |T6| . We get for N = 0,1 



Sjvx = —2ilm J du^2' 



(37) 



together with the self-consistent equation 



(j) r N (uj) = U(u) a+(u)5 N0 - a r (u)5 m - 

1 



— a r (u>) J did 
where II (a;) = l/[u — A — cr(u)} and 



L) - U? + IT] 



(38) 



dJ- 



to — uo' + irj 

is a self-energy which will play an important role in the following. Defining 

M") = ^4 - ^nn + ^n( w ) 



a r (uj) a r (u>) 



nM 



a(tu) 



U(oj) 



a r {uS) a r (uj) a(u) 
which are quantities independent of r due to Eq . (|38|) , we obtain the integral equation 



(39) 

(40) 
(41) 



[u-Aq- <t(w)] ip±{uj) = ± 



a ± (o;) 



alu) 



- duj' 



to — u' + irj 



Since Ima(u) = —na(uj), the solution of this integral equation is 

A- 



Imip±(uj) = ^7r^|n(cj)| 2 
A 



where 



A- 



A 



du\H(u) 



(42) 



(43) 



(44) 



and the real part of ip±(u}) can be obtained from ( f42|) and the Kramers-Kronig relation 
although it is not necessary for the following. From Eq.(^) we can calculate Im<j) r N (uj) with 
the help of (JD,|I) and obtain from Eq.(^) a nonperturbative expression for the transition 
matrix elements 



— E 



oo 



2m- 



X. 



A 



-En = 2m— — . 

A 



(45) 
(46) 
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Inserting these quantities in the kinetic equation (|35"D and solving it we obtain for the 
stationary probabilities 



P? = \- , Pf = A+, (47) 

where A + + A_ = 1 has been used which ensures the normalization P st + P[ l = 1. Further- 
more, both probabilities are strictly positive. The solution (|47|) is the final result for the 
density matrix and will be discussed in detail in section 5.2. 

In order to calculate the stationary current If we relate it to the correlation functions 
C > (u) and C < (uj) as indicated in Fig.^ which gives 



TSt 



-ie f du {a+H^H + a'^C^uj)} . (48) 



This relation is exact in the limit of a very large number of transverse channels. Otherwise 
the current will also depend on correlation functions involving more than two charge transfer 
operators e ±l< ^. 

The correlation functions can now be calculated from the diagrams shown in Fig.^8 
where we have used the same criterion as for the calculation of the density matrix with one 
exception. If a vertical line lies between the external vertices we allow for a cut through at 
most one solid line. Thereby we have used the fact that the external vertices also create one 
electron-hole pair excitation in the system. Formally we can say that a vertical line between 
the external vertices will in addition always cut the virtual line connecting the external 
vertices. The sum of all these diagrams gives 

<?» = mm \pfu{ u ) - y £ pp f dJ n(o,)l (49) 

[ JV=0,1 r=L/R LO LO +17] j 

C<(uj) = -2ilm I -PfU(coY + y y P S J [ du' U( 



M* (so) 



which yields with Eqs.(|^,[|l|),([|3|) and 



C > (lu) = -2ma-(Lu)\Tl(uj)\ 2 (51) 
C < {uj) = 2ma + (u)\Ii(uj)\ 2 . (52) 

Furthermore, we obtain for the spectral density (|i~8|) 

A(u) = a(u) |n(cu)| 2 (53) 

which is normalized to unity / duA(u) = 1. 

In summary we can express all our final results in terms of the spectral density 

^ = ^ ^j£M A(tf ) (54) 

Pf=f^?4P^) (55) 
J a{uj) 
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if = l^/^s ^y^ Mt^M - Ml (56) 
c>h = -arf g^£M ^H (57) 

a{uj) 

C<(u)=2 ^ a ^ )f ^ Uj) A(u) (58) 

where we have used the relation af(u) = a r {uj)f^{uj). 

These results satisfy conservation laws and sum rules. Current is conserved, i.e. J2r I? = 
0, and is zero in the equilibrium case when fi r = 0. All probabilities are positive, the 
spectral density is normalized to unity and we get the correct relationships between the 
correlation functions and the spectral density in the equilibrium case. The classical result 
can be recovered from ( |5"4] - |5"8| ) if we use the lowest order approximation in cto for the spectral 
density 

A(°\u) = 5(u- A ). (59) 

We conclude this section by the observation that quantum fluctuation effects due to energy 
renormalization and broadening manifest themselves in the spectral density via the real and 
imaginary part of the self-energy a(u) given in Eq. (p9|) . This quantity will be the main 
subject of the detailed discussion in the next section. 



V. RESULTS AND APPLICATIONS 



5.1 The spectral density 

Using Eq.(^) we can evaluate the self-energy u(u) given by Eq.(|39"D as 

a{uj) = — «q [R(bJ — H r ) + U(UJ — fi r )} 



(60) 



where 



R{uj) 



D(lu)uj 



^(^)+^(l 



E, 



c 



*2ttT' 



2ixT' 



2Reip{i 



M s 
2ixT' 



I(uj) = D(uj)-nuj coth 



l 2T' 



(61) 
(62) 



ip is the digamma function and we have chosen a Lorentzian cutoff function D(u) = E c / (u 2 + 
Eq) with a band width given by the charging energy Ec- For Ec ^> \oj — [i r \ T, the 
self-energy cr(u) behaves marginally 



aiu) 



E 



Qr 



2(u — fi r ) ln(- 



E. 



c 



\u — jl r 



(63) 



Since the imaginary part is linear in u, the broadening effects manifests themselves in a very 
unusual way. 
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We will now anlayse two limiting cases, either T \uo~ fi r \ <C Ec or \oj — fi r \ < T -C Ec- 
Furthermore, we assume for simplicity and with regard to the applications we are discussing 
in section 5.2 and 5.3 that all chemical potentials of the reservoirs are zero. 

In the first case, i.e. T <C \u\ <C E c , we can use Eq . (|63"D and obtain from Eqs.(^3|)and 
( |55D for the spectral density 

A(u>) - M A( ^ M , T « \u\ « E C (64) 

where A and a are the renormalized parameters 

AM= l + 2a?in(ff) ' 1 + (65) 

a(u) = ° -s - . (66) 

l + 2a ln(§) 

Approximately, A(u;) will have a maximum value at Ao = A(Ao) with a broadening of the 
order of 7rA a! where a = a(A ). For na <C 1 this broadening can be neglected and we 
obtain 

A ° l + 2a ln(g) ' a ° l + 2a ln(g) ^ 

for the renormalized gap and the renormalized dimensionless conductance. This result co- 
incides with the RG-analysis in ref. and for small «o where Aq has been replaced by Aq 



on the r.h.s. of Eq.(^) with the results of ref. ||11|| . This means that the leading logarithmic 
terms are included in our diagram series. 

For \u>\ < T <C Ec, we can approximate R{uj) by 2co>ln (J|jr) and obtain 

a, \ ~ A o coth(^) 

Aiuj) = —— ■ = — , \uj \< T <C En 68 

V ' A (^-A ) 2 + (7ra ^coth(^)) 2 ' y J 

where the renormalized parameters are now 

A = A _ = a 

° l + 2« ln(|^) ' a ° l + 2a o ln(i§0 1 ' 



which is the result (D) quoted in the introduction. If A < T, the spectral density 
has approximately a maximum at Ao with a broadening of the order of naoT. Again, for 
na <C 1 this broadening is small compared to A . However, the results ( j68| ) and fl6"9|) are 
independent of the value of A and can be used always when we need the spectral density 
only for |u;| < T. This is indeed the case for the calculation of the differential conductance 
in linear response as will be shown in section 5.3. 

As long as we are only interested in the maximum point A of the spectral density we 
can conclude that it is given by (|67|) for A ^> T and by (|69|) for A < T. This has also been 
proposed in ref. |T3f (in addition, we can see here that the correct replacement to get 
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from (|67| ) is given by Ao — > 2nT on the r.h.s. of (0))- Furthermore, by using our solution 
(|53|) for the spectral density, we can also describe the crossover regime between A < T and 
A ^> T and we are able to account for the finite broadening if ira approaches unity The 
latter case is of interest, since it is possible to achieve experimental parameters like 7^ ~ 3, 



a ~ 0.5 |42|] which gives 7ra ~ 0.75 by using Eq.(^). On the other hand, when ira <C 1, 
the effects of energy renormalization can become important when the temperature is very 
low. E.g. for 2aoln(^^) = 0.1, na <C 1 requires a ~ 0.01 according to Q6"5| ) which gives 
a temperature of the order of O.OOlEc- 

5.2 Charge fluctuations in the single electron box 

In the equilibrium case where /i r = 0, the SET-transistor becomes equivalent to the sin- 
gle electron box. The average excess particle number can be calculated from ([12]) and (|55|) 



N = J duf(uj)A(uj). (70) 

Within the classical approach given by Eq. (|59|) , one obtains 

N d = /(Ao) (71) 

where the bare gap Ao = Ec(l — 2C g V g ) can also be expressed in terms of the external 
gate voltage. Thus, N cl (V g ) shows a step at = (or A = 0) which is smeared by 
temperature. 



As can be seen from Fig.|19] and Fig.|20| there are clear deviations from the classical result 
if «o increases or if the temperature T is decreasing. To estimate this, we have neglected 
effects from the finite broadening, i.e. 7r«o "C 1, and have assumed that_ the energy of 
the ground state and the first excited state lie symmetric (with distance around the 



degeneracy point where V(l) = V(0) = [0. In this case the partition function reads 

Z= 2e-* cosh(^) (72) 
and we obtain for the average excess particle number N = n x — T-^^ \nZ (see Eq.(|7|)) 

From ( |BTD and flBTf ) it follows for Ao > T as well as for Aq < T that |^ = ^ where we 
have used «o < 1. Thus we find 



-(1 5 tanhl 

2 V l + 2a ln(|^) V 2T(1 + 2« In 



N=-(l- , — T-TE^ ^nhi — -—2——— )) , A <T , 7ra «l (74) 



N=l(l~ , , n _ \__, Er , ) , A »T , 7ra «l (75) 



2 V |A |l + 2a ln(g 



For A < T given by Eq. (|69|) , this means that we get significant deviations from the classical 



result © if 2a ln(ff) ~ 1. Thus by increasing a or decreasing T, quantum fluctuations 
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become more important. For Ao S> T, Eq.([75|) coincides with the result obtained in []11]] for 
the T = case. 

The slope at A = A = is given by 



dN 



1 



|A =0 



4T(l + 2a ln(|^))2 



(76) 



which differs from the classical result by the logarithmic term. 

Thus, we conclude that the presence of coherent higher order tunneling processes can be 
identified by an anomalous temperature behaviour of the slope of the Coulomb staircase at 
the degeneracy point. 



Fig.pi|a shows a comparison of the fits (|74| , |75|) with the correct result obtained from 
Eq.flTOD and (|53|). Both line shapes agree quite well for txq.q <C 1 where finite life-time 



broadening effects are not important. Fig.|21|b shows the same comparison for nao ~ 1 and 
we can see here clear differences for Ao > T whereas for Ao < T the approximation ( |T4|) 
still seems to be reasonable. 



5.3 Conductance oscillations in the SET-transistor 



In the linear response regime we obtain for the current from ([56 

-It = G(V L -V R ) 
with the conductance G given by 



TSt 

1 R 



(77) 



G = — -Art 2 / du HK ; v J A{u)f'{u) 
n J a{uj) 



(78) 



The derivative of the Fermi function restricts the integration variable to the regime u < T. 
Thus we can use (1681) and obtain 



(ko- 



uo/T 



a uj coth (^t) 



q = ^ 2?r 2 QqQq 

" h a J " sinh(cj/T) ( u - A ) 2 + (tt5 uj coth (^)) 2 



(79) 



where a = + aft and A are given by Eq.(^). 

Two analytic results can be obtained from this formula. First, the maximum conductance 
at A = A = is given by 



G, 



e 2 

T 2n ~ — 

n «q 



71 (tT« ) 2 -1. 

arctan 

2 v 2a rr 



(80) 



and secondly, the integral of the conductance over the gap A (or equivalently over e 2 ^V g ) 
is equal to 



dA G(A ) 



A a a 



a «o 



T. 



(81) 



The broadening 7 of the conductance peak then follows from 
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^ JdA G(A ) _ 7r 3 T(l + 2a In (ggQ) 
^moa; 7r — 2 arctan( 2 ^ ?r — j 



In the regime 7rao <C 1, (pOl) and (1521) reduce to 



^ ~ T ~W ' l + 2a ln(^) (83) 
7 = yT(l + 2a ln(^)) (84) 

These results lead to the following predictions. In the regime 2a hi(^0 <C 1, G max is 
constant and 7 is proportional to T. This is the classical result for sequential tunneling. For 
2«oln(^^) ~ 1 and na <C 1, G max and 7 contain terms logarithmic in the temperature, 
which are an indication for energy renormalization effects due to higher-order tunneling 
processes. For T — ► we have 

G max , 7~TlnT (85) 

i.e. the maximum value as well as the broadening go to zero. For 2aoln ( J^) ~ 1 and 
Ttao ~ 1, we have to account for energy renormalization effects as well as effects from finite 
life-times which lead to the complex formulas (|80| ) and (|82|). Fig.(|22|) shows the conductance 
versus A for several temperatures calculated numerically from formula (|78f) . The result 
demonstrates the predicted behaviour. 

For the conductance away from the degeneracy point Ao = 0, we can make the following 
analytic analysis. For Aq < T and nao < 1, we can replace the last fraction in Eq. (|79|) by 



5{uj — Ao) since the broadening is of the order ira^T <C T and the function sin ^^ T ^ varies 
on a scale of uo ~ T. Thus we obtain 

G = ^ gfei • A °/ T , A <T , ™ «1 (86) 

h & sinh(A /T) ' K ' 

which is the classical result but with renormalized parameters Ao and 5q. 

For Ao 3> T and 7r<So < 1 , we can replace the last denominator in Eq. ([79]) by 1/ Aq since 
w < T < A and iia u) < na^T < T < A . This gives 

G = J^o< ( £) 2 , A »T , 7ra <l (87) 

where the renormalization factor (1 + 2cto In has dropped out. Thus in the Coulomb 

blockade regime we recover the usual expression of inelastic electron cotunneling p6|Ji7[] 
provided that tcq:q < 1 and A ^> T. Therefore, the influence of resonant tunneling processes 
seems to be easier observable at resonance where A < T. Here we expect significant 
deviations from sequential tunneling already in the regime 2a In (7^) ~ 1 and effects from 
finite life-times (see Eq. (p0|)) for 7ia ~ 1. 

Finally, we will compare our result (|79|) with the one obtained in |20] which also coincides 
with other approaches |T^,|I^,^1|,0] . It is given by 



18 



r _ e\ 2 a o a o f , u/T a u coth(^) 

T h a J sinh(^/T) (cu- A ) 2 + (7ra AoCoth(|i)) 2 ' 1 ' 

Here the energy renormalization effects in «o and A have been neglected and the integration 
variable u is replaced by A in the broadening part of the last denominator. This corresponds 
to the introduction of a constant and finite life-time into the usual expressions of inelastic 



electron cotunnneling Hl6| , |T7|l which regularizes the integrals. This is only justified in the 
regime 2aoln(J^) 1 where the renormalization factor (l + 2a ln {■§^))~ 1 and the broad- 
ening 7r<So^ coth(^) ~ TtaoT <C T are unimportant. Thus the replacement uj — > A within 
the broadening doesn't change the result for sequential tunneling significantly (A < T) and 
for A ^> T Eq. (|88| ) leads to the result for electron cotunneling since the numerator of the 
last fraction of Eq. (fT9| ) has not been changed. However, for 2a ln(^^) ~ 1 or nao ~ 1, 
Eq.(|38"D can no longer be used at resonance as is demonstrated in Fig.( 



VI. CONCLUSIONS 

In this paper, we have aimed at presenting a detailed theory of quantum fluctuation 
effects in transport through small metallic islands with strong Coulomb interaction. With 
the help of a diagrammatic technique in real-time space, we have identified and evaluated 
the contribution of correlated higher-order tunneling events. Assuming a wide junction with 
many transverse channels, we have allowed different electrons to tunnel an arbitrary number 
of times coherently between the leads and the metallic island. Using the two charge state 
approximation we have included in a closed analytic form sequential tunneling, inelastic 
electron cotunneling and resonant tunneling processes. 

From a theoretical point of view it has turned out that the effects of quantum fluctua- 
tions can be understood very clearly by investigating the spectral density which describes 
the charge excitations of the system. In the classical regime, the spectral density has a 
sharp maximum at the gap energy Ao which is the difference of the Coulomb energies of two 
adjacent charge states. In the quantum regime, the maximum point renormalizes to Ao and 
we obtain a finite broadening which can be estimated to be of the order of t^oiq max{A , T} 
where o.q is the renormalized dimensionless conductance of a single barrier. Both features 
together with the complete form of the spectral density is described by the real and imag- 
inary part of the self-energy a(u). It contains an anomalous dependence on energy and 
temperature which leads to a variety of unexpected features in the line shapes of several 
experimental quantities. 

To estimate the experimental consequences of quantum fluctuations we have calculated 
the average charge of the single electron box and the linear conductance of the SET-transistor 
as function of the gap energy Aq or equivalently the external gate voltage. For the average 



charge we have compared our results to previous investigations where renormalization 



group techniques have been used to study the zero temperature case. Using the temperature 
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as a cutoff one can calculate the average charge from these theories in the two limiting case 
A > T or A < T. These results agree with our solution in the case 7ra ^ 1, i-e. in 
the regime where finite life-time effects are not important. For 7ca ~ 1 there are significant 
deviations at least for A > T. Furthermore our complete solution is capable of describing 
the complete crossover from A < T to A 3> T and shows that the temperature has to be 
introduced into the renormalization of the system parameters by the replacement Ao — > 2itT. 

For the conductance in the linear response regime we have seen that the classical descrip- 
tion of sequential tunneling near the resonance (Aq < T) is only valid for 2a rn(J^) < 1 
and na = na /(l + 2a ln(^^)) <C 1. Nowadays experiments are not restricted to this 
regime. Therefore we expect that resonant tunneling of coherent higher order tunneling pro- 
cesses should be observable by a measurement of the line shape of the conductance peaks 
as function of «o an d temperature. The renormalization of the system parameters A and 
«o is important for 2a ln(^;) ~ 1 and leads to an anomalous logarithmic temperature 
dependence of the conductance peak and the broadening. For T — > the conductance maxi- 
mum decrease like 1/lnT and the broadening increases proportional to Tin T. Furthermore, 
for nao ~ 1, the influence of finite life-times becomes very important and leads to a very 
significant flattening of the Coulomb oscillations. Both effects are important and we have 
seen that it is very difficult to separate them at realistic temperatures. 

Finally our approach describes also the conductance in the Coulomb-blockade regime 
where transport is dominated by inelastic electron cotunneling. Our analytic formulas give 
the correct crossover from resonant tunneling at the degeneracy point to inelastic cotun- 
neling. We have seen that cotunneling persists in the regime A ^> T and 7ra < 1. For 
nao ^> 1 cotunneling as well as the validity of our approach will break down. In this regime 
the charge will no longer be well defined and one should go beyond the two charge state 
approximation. 
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FIGURES 



FIG. 1. Equivalent circuit for the SET-transistor. 

FIG. 2. Example of a diagram in time space for the probability Pjv(i)- Time is increasing from 
left to right. The dots represent the vertices which are connected by horizontal lines (propagators), 
solid lines (leads) or wiggly lines (island). €j,Ej are the corresponding energies, rj are the lead 
indices and N denotes the excess particle number. 

FIG. 3. Graphical representation of (a) the current I r through lead r and (b) the correlation 
functions C > , C < in time space. The dots indicate the external vertices. Further internal vertices 
together with their connections are not indicated. 

FIG. 4. Graphical determination of the energy denominators. The energies of all lines con- 
tribute to AE which are cut by an auxiliary vertical line. 

FIG. 5. If an auxiliary vertical line cuts the virtual line connecting the external vertices, the 
energy u contributes to the energy denominators. 

FIG. 6. Graphical representation of the correlation functions in energy space. 

FIG. 7. Self-consistent equation for the stationary probability Pj$. The self-energy £ denotes 
the sum of all irreducible diagrams which can not be cut into two parts by an arbitrary vertical 
line. 

FIG. 8. Examples for irreducible diagrams. An arbitrary vertical line will always cut through 
some reservoir line. 

FIG. 9. First order diagrams for £ which correspond to the classical transition rates. 

FIG. 10. Replacement of a loop consisting of one lead and one island line by a single solid 
line. The arrow direction is the same as of the lead line. 

FIG. 11. Second order diagram for £oo which contributes to inelastic electron cotunneling 
when the two leads r and r' are different. 

FIG. 12. Higher order diagrams which contribute to resonant tunneling. 
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FIG. 13. Example of a diagram which is taken into account. All vertical lines cut at most 
through two solid lines. The cut through the two horizontal lines is trivial and is not counted. 

FIG. 14. Graphical representation of the self-energy S^vi within our approximation. 

FIG. 15. Self-consistent equation for <p r N {uj). 
FIG. 16. Definition of II. To each energy denominator one has to add u> in order to obtain 

n( w ). 

FIG. 17. Graphical representation of the relation between the current and the correlation 
functions. Here the line connecting the external vertices is a real one. All the other internal 
vertices are not indicated. 

FIG. 18. Graphical representation of (a) C > (uj) and (b) C < (uj) within our approximation. 

FIG. 19. Average charge as function of the gap energy for different values of ao- 
Ec = 1, T = 0.01 and (a) ao = , (b) a = 0.01 and (c) a = 0.1. Curve (a) is the Fermi 
distribution function which corresponds to the classical result. For increasing ao the deviations 
become more significant. 

FIG. 20. Average charge as function of the gap energy at zero temperature. At finite ao = 0.03 
there are clear deviations from a pure step function. 

FIG. 21. Average charge as function of the gap energy using (1) the correct result Eq.(67), 
(2) the fit (71) and (3) the fit (72). E C = 1,T = 0.01 and (a) a = 0.02 and (b) a = 0.2. When 
finite life-time broadening effects set on in (b) the deviations grow. 

FIG. 22. Conductance as function of the gap energy for various temperatures. 

E c = 1 , a = 0.02 and (1) T = 0.05 , (2) T = 0.02 , (3) T = 0.005 and (4) T = 0.0005. 

FIG. 23. Conductance as function of the gap energy using (1) the correct result (76), (2) the 
fit (86) and (3) the fit (84). E c = 1 and (a) T = 0.001 , a = 0.02 and (b) T = 0.05 , a = 0.6. In 
both cases are clear differences to the classical result (86). (a) can be described by renormalized 
parameters given by (84) whereas in (b) finite life-time broadening effects dominate. 
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